% Specify parameters

param.L = 0.2; % labor endowment
param.zeta = 4 ;% elasticity of substitution between domestic and imported final goods
param.vartheta = 0.4; % elasticity of substitution between tradables and non-tradables
param.sg = 3.8357;% demand elasticity
param.eps = 2.3777;% the elasticity of substitution between domestic and foreign varieties in production function
param.gamma = 0.6079;% output elasticity of materials
param.eta = 0.382;% eta, sensitivity of foreign input price-index to mass of countries sourced

% For sensitivity analysis
if eta35glob==1 || eta35NGF==1
    param.eta = 0.35;
end
if eta045glob==1 || eta045NGF==1
    param.eta = 0.446;
end
if eps27glob==1 || eps27NGF==1
    param.eps = 2.7;
end
if eps22glob==1 || eps22NGF==1
    param.eps = 2.2;
end

sg = param.sg;  g = param.gamma; eps = param.eps; eta = param.eta;

% Normalizations
param.z  = 0.001; % level of price index of foreign input bundle
param.qD = 1; qD = param.qD; % quality of domestic input bundle
param.AN = 0.05; % productivity in non-tradable
param.aN = 0.35; %location parameter for N in utility
param.pM = 1; % price of imported manufacturing final good
param.tau = 0; param.bl  =1; % export iceberg cost and scale parameter of distribution of foreign demand across countries 
knorm = sg^(-sg)* (sg-1)^(sg-1) * (1-g)^((sg-1)*(1-g)) * (g * param.qD)^((sg-1)*g); % see calibration appendix
param.knorm=knorm;

% Firm heterogeneity
if glob==1
    disp('Baseline model')
    load parameters_baseline.mat
elseif nafta==1
    disp('NAFTA counterfactual')
    load parameters_baseline.mat
elseif ngf==1
    disp('No global firms model')
    load parameters_NGF.mat  
elseif monly==1
    disp('Model of importing')
    load parameters_importing.mat
elseif eta35glob==1
    disp('Baseline with eta=0.35')
    load parameters_eta035.mat 
elseif eta35NGF==1
     disp('No global firms model with eta=0.35')
     load parameters_ngf_eta035.mat
     x(5) = 1.0000e+40;
elseif eta045glob==1
    disp('Baseline with eta=0.45')   
    load parameters_eta045.mat
elseif eta045NGF==1
    disp('No global firms model with eta=0.45')   
    load parameters_ngf_eta045.mat   
elseif eps27glob==1 
    disp('Baseline with epsilon=2.7')   
    load parameters_eps27.mat
elseif eps22glob==1 
    disp('Baseline with epsilon=2.2')   
    load parameters_eps22.mat
elseif eps27NGF==1
    disp('No global firms model with epsilon= 2.7')   
    load parameters_ngf_eps27.mat
elseif eps22NGF==1
    disp('No global firms model with epsilon=2.2')   
    load parameters_ngf_eps22.mat
end

param.std_phi = x(6); % dispersion of efficiency
param.a1 = x(1); % controls average fixed cost of importing
param.a2 = x(9); % controls correlation between fixed cost of importing and efficiency
param.a3 = x(7); % controls dispersion of fixed cost of importing
param.a4 = x(2); % controls average fixed cost of exporting
param.a5 = x(10); % controls correlation between fixed cost of exporting and efficiency
param.a6 = x(11); % controls correlation between fixed cost of importing and fixed cost of exporting
param.a7 = x(8);  % controls dispersion of fixed cost of exporting
FM = x(3); param.FM = FM; % fixed cost of being an importer-only
FX = x(4); param.FX = FX; % fixed cost of being an exporter-only
FXM = x(5); param.FXM = FXM; % fixed cost of being an importer-exporter
param.FM = FM; param.FX = FX; param.FXM = FXM;
param.theta = x(12); theta = param.theta;

param.beta =0.5; beta = param.beta; % beta is location parameter in the CES aggregator of materials; beta is the weight on the domestic variety
betastuff = beta ^ (eps / (eps-1) * g * (sg-1)); betastuff2 = beta ^ (eps / (eps-1) * theta * g * (sg-1)); betastuff3 = (beta / (1-beta))^(1/eta * eps/(eps-1)); % constants used to compute variable profits

if glob==1 || eps27NGF==1 || eps22glob==1  || eps22NGF==1  || eps27glob==1 || eta045glob==1  || eta35glob==1 || eta35NGF==1   || nafta==1  || eta045NGF==1
    param.number_firms = 200000; 
    gridsizez = 200;  gridsizey = 800; gridsizes = 800;
    mingridz = -3.7; maxgridz = 4.5; mingridy = -40; maxgridy = 100;
    mingrids = log10(1e-08); maxgrids = log10(1-1e-08);
elseif  ngf==1 
    param.number_firms = 200000; 
    gridsizez = 200;  gridsizey = 800; gridsizes = 800; 
    mingridz = -3.7; maxgridz = 4.5; mingridy = -40; maxgridy = 100;
    mingrids = log10(1e-08); maxgrids = log10(1-1e-08);
elseif monly==1
    param.beta =1; beta = param.beta; betastuff = 1; betastuff2 = 1; betastuff3 = 1;
    param.number_firms = 49003; 
    gridsizez = 1200;gridsizey = 100; gridsizes = 1100;    
    mingridz = -4; maxgridz = 30; mingridy = -4; maxgridy = 50;
    mingrids = log10(1e-08); maxgrids = log10(1-1e-08);
    param.aN = 0.1; param.AN = 0.0175;    
end

if eta35glob==1 || eta35NGF==1 || eps27NGF==1
    maxgridz = 6; 
end

param.betastuff = betastuff; param.betastuff2 = betastuff2; param.betastuff3 = betastuff3;
number_firms = param.number_firms;
maxshare = 1-1e-08; minshare = 1e-08; % Bounds for interpolated shares

Opt_options = optimset('Display', 'iter'); % Optimiation options for devaluation function
Opt_options.MaxFunctionEvaluations = 500;

param.dz = 1; 
param.dAN =1;
param.dbl = 1;